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A scheme for measuring complex temperature partition functions of Ising models is introduced. In 
the context of ordered qubit registers this scheme finds a natural translation in terms of global operations, 
and single particle measurements on the edge of the array. Two applications of this scheme are presented. 
First, through appropriate Wick rotations, those amplitudes can be analytically continued to yield esti- 
mates for partition functions of Ising models. Bounds on the estimation error, valid with high confidence, 
are provided and shown to be compatible with previous results. Interestingly, the kind of state prepa- 
rations and measurements involved in this application can in principle be made "instantaneous", i.e. 
independent of the system size or the parameters being simulated. Second, the scheme allows to accu- 
rately estimate some non-trivial invariants of links. A third result concerns the computational power of 
estimations of partition functions for real temperature classical ferromagnetic Ising models on a square 
lattice. We provide conditions under which estimating such partition functions allow to reconstruct scat- 
tering amplitudes of quantum circuits making the problem BQP-hard. Using this mapping, we show 
that wavefunction overlaps for ground states of quantum Hamiltonians, which can serve as a witness to 
quantum phase transitions, can be estimated from classical Ising model partition functions. Finally, we 
show that the ability to accurately measure corner magnetisations on thermal states of two-dimensional 
Ising models with magnetic field leads to fully polynomial random approximation schemes (FPRAS) for 
the partition function. Each of these results corresponds to a section of the text that can be essentially 
read independently. 

1 Introduction 

Statistical Mechanics provides formal recipes to study interacting many-body systems. Quantities that can 
be experimentally probed, such as the free energy or the specific heat, can in principle be derived straightly. 
More often than not, though, computing these quantities turns out to be impossible. As can be seen from 
very idealised systems, our ability to actually apply these recipes is very limited. During the last ten years, 
significant efforts have been devoted to investigating whether quantum mechanics could help in this respect. 
Various methods, all involving the superposition principle, have been proposed to compute the Jones poly- 
nomial at particular values of its variable [IJ, partition functions of classical statistical models lEllS, the 
Tutte polynomial L4J, or more generally to contract tensor networks [3 |. 

In this work, we will mainly focus on a collection of classical two-level systems, each attached to a 
fixed position corresponding to a vertex of some lattice A, with edges E{A). The state of a particle located 
at vertex i is associated with a number taking value in { — 1, +1}. The energy of the system is given by 
an Ising Hamiltonian function, associating an energy with each classical configuration of the system a\: 
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The first sum in this equation runs over all vertices of A. The quantity hi models represents some local 
field felt by a spin located at position i. The second sum represents interactions between pairs of neighbour 
particles (edges of the lattice). The strength and sign of these interactions may vary from pair to pair. 
This model was introduced by Lenz as an idealisation of systems where magnetic interactions prevail [3. 
Although innocent looking, it exhibits an extremely rich structure. On a regular lattice, close to a phase 
transition, its long range behaviour is similar to that of very interesting field theories |6| while the problem 
of computing its partition function, 

Z(/3) = 5^exp[-/3i/(M)], (2) 
W} 

belongs the NP-hard complexity class |7l. 

It is the purpose of this paper to present schemes that allows to accurately estimate Z{f3) for imaginary 
values of /3 (Section [2]), through manipulation of a suitable quantum mechanical system. Quantum circuits 
for this task have been previously proposed in Ref.|8|. However with our scheme, we will see how to eval- 
uate partition functions of real systems, through analytic continuation (Section]?]). As we shall see, the kind 
of preparation and measurement necessary for this estimation can in principle be made in constant time, i.e. 
independent of the system size or the parameters being simulated. This feature is particularly appealing in 
view of possible practical implementations. We will then see that imaginary temperature partition functions 
are interesting in their own right, because they provide non-trivial invariant of knots (Section |5]). Section 



6.1 deals with computational complexity issues. We investigate the (quantum) computational power of the 



Ising model, and show how the ability to estimate real temperature partition functions of this model allows 
to efficiently simulate a quantum computer. One application of this is the estimation of the wavefunction 
overlap between ground states of a quantum Hamiltonian in the vicinity of a quantum phase transition. We 
also show that some much simpler tasks have computational power. In particular, the ability to detect cor- 
ner magnetisations of disordered Ising models leads to fully polynomial random approximation schemes 
thereof. Many of the quantum algorithms presented here involve repetitions of either constant depth or 
linear depth circuits and moreover many of the operations can be performed without individual qubit ad- 
dressability. This is potentially a real boon to experimental implementations in architectures such as trapped 
atoms in optical lattices or superconducting qubit arrays where individual addressing is not so easy but many 
qubits are available. In additional some of the circuits provide for a trade off in space and time, i.e. one can 
perform either constant depth circuits in d + 1 spatial dimensions or linear depth circuits in d dimensions. 
Constant depth quantum circuits have attracted attention since the discovery of simple examples (depth- 1 
circuits) that are expected to be difficult to simulate classically |[9l- Furthermore, there is some evidence 
that fault tolerance thresholds could be improved for constant depth (or more generally logarithmic depth) 
quantum circuits |[T0l[TT1l . 



2 Complex temperature partition functions 

We wish to study a classical system defined on some d-dimensional lattice A. For that purpose, we consider 
an associated situation, where a two-level system is located on each vertex of A. The computational basis 
for each quantum particle, {|+), |— )}, will be associated with classical individual spin configurations. Our 
construction relies on controlled phase gates acting on nearest neighbours, that is, elements (A:, /) of E{A), 
the set of edges of the lattice. Their action is best described in computational basis: 

Ck,i:\ak,a'i)^e"t"'-^^'''"''i^ak,a'i). (3) 
Importantly, these phase gates all commute with each other: 

\/{k,l),{x,y) eE{A), [Ck,i,C,,y]^0. (4) 
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Obviously, each function 0(/e,^> can be expressed as 



s'(t[- 

s=±ls'=±l 

With the definitions A^/e = \ Y.sk,si 4^k,i{sk, si), Jk,i = lJ2sk,si 4^k,i{sk, si)skSi,hk = iJ2sk,si 4^k,i{sk, si){sk+ 
5/), we see that a collective action of controlled phase gates across all edges of the lattice can be described, 
in computational basis, a^ 

n Cj^'^i]^\ak) = oxpliaYf^k + ic^Y^^^^^^^ X] Jk,l(^k(^l] YlWk)- (5) 

{k,l)eE keA keA keA {k,l)eE keA 

In particular, if each quantum particle is initialised in the state 

|+,)^i=(|+) + |-)), (6) 

we see that the mean value of a product of phase gate operators takes the form of a partition function at 
imaginary temperature ia: 

{k,l)eE {a} 

with H of the form given by Eq. ([!]). 

It is actually possible to get partition functions of a classical {d + 1) -dimensional system through evolu- 
tion of a d-dimensional quantum system. For that, we use two additional kinds of gates besides the controlled 
phase gate. The first kind is single qubit rotations: 

Uk'-\+) cos6>fc|+) + sin6>fc|-), 

Uk'-\-) -sm9k\+) +cos9k\-). 

(8) 

As discussed in Appendix [A] other choices are possible. The second is single qubit phase gate: 

Pk{^k)-Wk)^e'''''''Wk)- (9) 

Next, we observe that the matrix elements of Uk can be expressed in exponential form for almost all values 
of the parameters 9]^: 

{(^k\UkWk) = exp [Jtdkak + - i'^^k + B(9k)] , 0^ ^ : G Z}, (10) 

4 = -^ln(tan0,)-z^, (11) 
and 

ln(cos(gfe)) ln(sin(gfe)) tt 
B{Ok) = ^ + ^ + . (12) 

These individual rotations {Uk, /c G A} are applied on all lattice sites simultaneously. For bookkeeping, it 
is convenient to assume there is an external clock recording the moment t where simultaneous rotations are 



where 



^Note to a reader interested in reproducing the calculations: the identity 6cra' = ^~^f^ been repeatedly used. 
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applied, and ticking at exactly this time. There is nothing particular to this clock, it is just a way to label the 
change of variables necessary to describe the action of the Uk gates: 

l[Ukit)\ait)) ^ G{t) exp[^4(tK(f)<7fc(t+l)+i^^(<7fc(t+l)-afc(f))]|a(t+l)), (13) 

keA {cr(t+l)} keA keA 

where G{t) = e^EkeA ^(^^ W))- 

Now let us consider a d-dimensional lattice A of particles each prepared in the state ([6]). Let us assume 
that a layer evolution operator 

^w=n^^(-i)nf^^w n c'.^^wn^fcij)- ^4) 

keA keA {k,l)eE keA 



is applied (m — 1) times on this initial state, leading to the final state nj^i^ ~ (see Fig{l]). 

The overlap of this state with the initial state |+x '^') takes again the form of an Ising partition function, 
but now defined on an enlarged lattice A = Ax{l,...,m}: 

m—1 

A{a,e) ^ n ^(m-t)|+?l^l) ^-^eM-Hia)], (15) 



t=l 



where, 6 denotes collectively all individual rotations performed on the system, and where, up to an additive 
constant J2t^i 1^ G{t), the classical hamiltonian H is given by 



m—1 



t=i keA t=i {k,l)eE t=i keA 



Eq.(15) is proven by inserting identity operators and identifying single-particle quantum basis states |zt) 
with single particle classical spin configurations \a): 



m—1 m—1 



2^ 

t=l a(l)...cr(m) t=l 



3 Implementation 

At the core of the discussion held in the previous section lies the ability to measure the scalar product 
between n-particle states |$) and |^^). We will describe two measurement protocol addressing this problem. 
The first is the simpler and allows to detect | ($|^^) p, while the second truly yields ($|^^). 
Protocol 1 

1. Prepare an n-particle system A in the state |^^), and an n-particle system B in the state |$). 

2. Prepare an ancillary register Rofn qubits in the state \GHZ) = -^(| + ...+) + | — . — )). 

3. Perform a bit- wise controlled swap gate with each qubit Rj of the register as a control and Aj, Bj as 
targets, i.e. if qubit Rj is in the state |— ) then apply SV^AP (Aj ^ Bj). We get 

^ (l + . . . + I - . . . -)r\^)Mb) 
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Figure 1: Example of the quantum algorithm on a ID chain of qubits to compute the partition function 
of a 2D classical Ising model at imaginary temperature, (a) The quantum algorithm begins with qubits in 
the chain initialized in state 1+^;) and proceeds with alternating sequences of parallel nearest neighbour two 
qubit gates Ck (t) diagonal in the computational basis { | ±) } and parallel local rotations Uk (t) (supplemented 
by single qubit phase gates), (b) The corresponding classical Ising model with spatially dependent horizontal 
and vertical bond strengths and local magnetic fields. 

4. Measure the first n — 1 qubits of R in the basis = ^'^^^^'^ ^ Denote rrij = ±1 the 
(equiprobable) outcomes of measurement on register qubit j and define x = ^j- ^^e state for 
the last qubit of the register and the system AB is 

+ (17) 

5. Measure the Pauli operator of the last ancillary qubit i?^. The expectation value is 

«) = (-lf|($|v|>)|2. (18) 

Protocol 2 

1. Prepare an n-particle system A in the state |$). 

2. Prepare an ancillary register Rofn qubits in the GHZ state 

3. Evolve the qubits in register A conditioned on the state of the ancilla to prepare 

This can be done by replacing all instances of quantum gates in the evolution of |$) |+f ^) |^) 
into bitwise controlled gate operations. The single qubit phase gates become controlled phase gates: 
I (+ 1 ® 1 Afc + I (- 1 ® • Similarly, the single qubit rotations become: | +) i^^ (+ 1 ® 1 + 
\—)Rk{~\ ® ^Afc- The collisional gates are controlled by one of neighboring ancillary qubits, e.g.: 
\+)Rx (+1 ® '^Aa:,Ay + \-)Rk {~\ ® CA^.Ay Such three qubit diagonal gates can be decomposed into 
at most 6 nearest neighbor controlled phase gates l[T2ll . 
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4. Measure the first n — 1 qubits of R in the basis { | ) } . Denote mj = ± 1 the outcome of measurement 
on register qubit j and let again x = ^]=i ^^^^^ ^^^^ qubit of the register and the 
system A is 

-^(|+,$) + (-iP|-,v|>)). (19) 

5. Measure cr^ on the last ancillary qubit i?^. The expectation value is 

K) = (-im($|*)]. (20) 



6. Repeat steps 1-4 but on the last qubit Rn measure instead the Pauli operator where the basis 
{|±^)} = {-^(1+) ±z|— ))}. The expectation value is 

«) = (-l)>^55[($|v|>)]. (21) 



We note that it is actually not necessary to prepare size n ancillary registers in a GHZ state for either 
measurement protocol, since one ancillary qubit making controlled swaps or controlled interactions like a 
serial tape head over the quantum registers would suffice. The penalty is a potentially linear slowdown and 
the need to transport the ancilla qubit over the register for every gate in the circuit. The \GHZ) state can be 
prepared in one plane using global, i.e. spatially homogeneous, pulses in the plane ||l3l[l4l[T5l. Futhermore, 
by coupling the quantum register with a common bosonic mode, \GHZ) states can be prepared in constant 
time |[T6ll . The idea is to place all the spins inside a high Q cavity (with decay rate with a resonance field 
frequency close to the transition between the qubit states and some other excited state. When the coupling 
between the field and qubits is spin dependent and dispersive (e.g. a differential light shift induced by 
polarization section rules or by spin dependent detuning) then the interaction is modelled as: 

j 

where Qz is the dispersive coupling strength. Then \GHZ) can be produced either using strong coupling 
with a quantised state of light or via a geometric phase gate using coherent state displacements. We outline 
the latter as follows: 

• Initialize all the spins in 1+^;) and the cavity mode in the vacuum state |a = 0). 

• Perform the following nine step interaction sequence: 

j j 

where D{a) = ^ coherent state displacement, and e~^^^^ is the unitary evolution gener- 

ated by Vz. When the parameters satisfy: gzT = 7r/2, and |a/3|(e~^^'^/^ + e"^'^/^) = 7r/4, then the 
cavity returns to the vacuum and the global rotation U = e~^^ ^| applied to the qubits. 

• Apply the global operation Hj ^ to the spins. 

The state of the qubits is then ^(| + + ...+) — i\ ) which is locally equivalent to \GHZ) and 

functions just as well for the simulation protocols above. The overall process fidelity, which measures how 
close the lossy process is to the target unitary U = 4 ^ satisfies L17il 

i^pro > 1 + ^ . 

2|5^| V 215^1/ 
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Note that this is a constant depth circuit thanks to the non-local coupling of the field to the qubits. Of course 
as the number of spins increases the size of the cavity must also increase, and the strength of the field, spin 
coupling decreases as 1/ VVol where Vol is the cavity volume. Consequently, there is ultimately a process 
time which scales as y/n where n is number of qubits. However, in practice this could be quite fast compared 
to a sequential circuit for generating \GHZ). 

Since the measurement of (^^|^^) is informationally more complete than that of the reader 

might wonder why we have bothered describe a separate procedure to measure the latter quantity. The 
reason is that the second measurement seems to be experimentally more demanding than the first. The 
information provided by the first kind of measurement is however certainly valuable and for most of the 
discussion to follow we we assume data is obtained from this kind of measurement. 

So far, we have considered the case of planar boundary conditions. If the classical system is periodic in 
space (i.e. the lattice A is periodic) then the above quantum algorithm is simply modified in the couplings 
Jk,i{t) to account for this. If the classical system is periodic in the time direction, then a few modifications 
are needed. To relate the measurement of the quantum system to the classical partition function, the bound- 
aries states \ak{m)) and |crjt(l)) must be identified. So rather than computing the scattering matrix element 
(+x '^'), where the unitary W is defined as = 11^^^^ C{m — t), as we have described so far, we 
want the trace: Tr[VF]. This is found by using the measurement protocol 1 but with the register A prepared 
in the completely mixed state The polarization measurements of the last ancilla of the register then yield 

the real and imaginary parts of ^ . Also note by the cyclic property of the trace, the phase gates Pk are 
no longer needed in the quantum evolution. 

Consider the implementation of this measurement for a 3D classical Ising model using a quantum reg- 
ister encoded in a plane. For protocol 1 three parallel planes are needed, one (the top plane) prepared in a 
\GHZ) state, and the centre (c) and bottom (b) planes both prepared in |+f The centre plane is prepared 
in ]l{k,i)^E C'^^/l+x'^') or evolved in Y\TJi >C(m - t)\+x^^^), and the subsequent C-SWAP gates between 
registers can be implemented in parallel bitwise between pairs (c/e, hk) using a sequence of at most 12 near- 
est neighbor collisional gates [12]. Finally the measurement of the top register only requires collecting the 
parity of measurement outcomes of n — 1 qubits in the bulk (without addressability) and an addressable 
measurement of Xn for one qubit on a corner. For protocol 2 two registers are needed: the top one prepared 
in \GHZ) state and the bottom prepared in l+f*^). During the quantum evolution all gates acting on the 
bottom register (say qubit hk) are to be controlled by the neighbouring qubit on the top plane (qubit tk). For 
a rotation gates Uk{t) this means to instead apply the controlled gate ® l^^ + \~)tk{~ \ ® Uk{t). 

Such a gate can be done using at most 3 controlled collision gates between and b^. For the two qubit 
gates Ck^i{t) we need to apply 1+)^^ (+| ® ^h^M + \~)tk {~\ ® Ck^i{t). This three qubit diagonal gate can be 
realized using using at most 12 collisional gates between nearest neighbors t^, hk and hk^hi. Since not all the 
gates now commute it is necessary to do this in two stages over non overlapping pairs of nearest neighbors 
in the bottom register. Measurement of the top register proceeds as for protocol 1. Regarding addressability, 
it is necessary to be able to address the different planes along z but addressability can be relaxed in the x — ^ 
direction. 

4 Partition functions 

The schemes of Section [2] can be used to provide estimates for real temperature partition functions of clas- 
sical models. We proceed by analytic continuation of the quantum amplitudes (or their modules) provided 
by the protocols described in Section [2j The general idea is to view the partition function as a polynomial 
of order linear in the system size whose coefficients are the same as the those obtained from the quantum 
amplitude estimation but with real instead of complex variables, and then to Wick rotate these variables. 
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Let a and 6 denote two complex variables, and consider a function F of the form 



Ni N2 

F:CxC^C:{a,9)^ F{a, 9) = Yl Yl c^l,^2e'^'''e^'''^ (23) 

Ui= — Ni U2 = — N2 

where A^i, A^2 < oc. Clearly, F is an analytic function, so the coefficients {ciy^^iy^} define F on the whole 
complex plane. If F is known for aj^ = a^-^'i) = 27r ji = . . . 2Ni,9j^ = 6>(-^2) = 27r ^, j2 = 
. . . 2A^2, then a Fourier transform yields 

Plugging this expression in Eq.([23]), one finds sums of geometric series. Summing them yields 

2Ni 2N2 

F{a,9) = E E - c^nW'Ho - 9j,), (25) 

jl=0j2=0 

where 

1 dn((2Ar + l)i)^ 

^ ^ 27V + 1 sm f 

Now consider the quantum amplitudes introduced in Section [2| in the case where hk{t)^ Jk^i{t) G 
{-1, +1}, Vfc G A, V(A:, /) G E{K), Vt = 1 . . . m, and where all 'vertical' couplings are set to -1. 
(For the case of non-uniform vertical couplings, see Appendix |A|) In that case, these quantum amplitudes 




are certainly of the form ( |23| ), with A^i , A^2 growing at most polynomially with the number of vertices of the 
classical model being under consideration. For suitable complex values of a, the probability amplitude 
A{a^9) of the d-dimensional quantum system can be put in correspondence with the real partition function 
of the {d + 1) -dimensional classical system. Namely, for 

1 ZTT 1 

= = 7 1^ \/ \ i l2^ji ' 9m = 2 1^ sin 2^ + ^ - - In 2, (27) 

one finds that A(a^,6>^) = el^l^^(^*) Z^^'^'^^li) /2\^\ . In the definition of a*, we recognise the familiar Wick 
rotation. The role of the other parameter, 9~^, is to analytically continue the unitary quantum mechanical 
transfer matrix, between successive times, to the (non-unitary) statistical mechanical transfer matrix. In 
summary, in order to get information about the partition function of a d-dimensional classical system, we 
estimate the probability amplitude A{a^ 9) for well-chosen values of a and 9. From the collected data, we 
reconstruct the dependence of the function A on its variables (q;,6>), as just explained. Finally, analytic 
continuation of the variables (a, 9) to the suitable values ( [27] ) yields an estimate for the desired partition 
function. 

Let us analyse the errors appearing when the values A{aj-^ , 9j^) are not known exactly but estimated by 



some quantities ^j-^j2' The identity (25) allows to get a priori error estimate. To simplify the discussion, 
let us start with the case where partition functions are estimated using a one-time-step protocol. Then, 
m = 0, A^i = poly(|A|) = N and A^2 0. Defining dLp"" max{|(/p^- - A{aj)\,j = . . . 2A^}, we see, 
through error propagation, that the error at inverse temperature /3, /S.A{il3) satisfies 

2N 

AA(z/3) < ^ \w^^\i^ - aj)\5i^\ (28) 
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In the limit of large values of /3, the r.h.s of this equation essentially behaves as S(f~^e^^, indicating that the 
measurement accuracy should grow exponentially small, with the inverse temperature and the size of the 
system, in order to maintain the error over our estimate for partition function below some fixed prescribed 
threshold. 

A bound on the error independent of /3 can also be derived easily. Indeed, for the Hamiltonians we are 
considering, the partition function can be written as 

N 

k/3 



k=-N 



where all coefficients are non-negative integers whose magnitude is at most 2^l^l (number of classical 
configurations associated with the system). It would therefore be sufficient to be able to estimate these 
coefficients with a relative accuracy of 2~^l^l in order to be able to reconstruct Z^^^^^{^) perfectly. The 



bound appearing on the r.h.s of ( |28) ) is independent of the actual values for the link couplings and magnetic 
fields of the precise Ising model being simulated. We therefore expect it to be pretty loose. 

To get a sharper understanding of how errors behave, we made some numerical simulations. In Figj2] 
we show how the error behaves by studying different quantities such as the logarithm of the partition func- 
tion, the energy and the specific heat. In particular we simulated a model with uniform couplings and zero 
magnetic fields and a model with ±1 couplings (with 50% probability) and uniform magnetic field. One 
can appreciate how, in the uniform case, the error over the partition function goes to zero for zero and in- 
finite temperature. In Appendix [A] we show that error over each Fourier coefficient is well behaved for 
large and small values of k (close to =bA^), but blows up for intermediate values k (close to 0). This fact is 
consistent with our numerical observations and the well known duality present in this model ifTSl . For the 
non-homogeneous case, we have found that the errors in the partition function starts by growing exponen- 
tially with /3, then remains constant. This observation is consistent with the fact that there is no known low 
temperature/high temperature duality relation. The errors we have found are also much larger. Our numerics 
indicate that, in the non-homogeneous case, the magnitude of the partition function is dominated by those 
corresponding to intermediate values of k, much more so than for the homogeneous case. 

Previous attempts at using quantum mechanics to compute approximations of partition functions exhibit 
errors comparable to ours. A quantum algorithm based on Fourier sampling was introduced in [19| to 
estimate partition functions and free energies of quantum Hamiltonians, which includes the classical Ising 
model in the case of all diagonal interactions. There it was found that the number Fourier components needed 
to be sampled scales polynomially with the lattice size, but in order to obtain a mulitplicative approximation 
of the partition function, the requisite accuracy of estimation of each coefficient scaled exponentially with 
the system size. An algorithm, based on using a quantum computer to contract tensor networks yields similar 
approximation scales |3 ]. Even preparing a quantum state which coherently encodes a classical thermal state 
of an Ising appears to be difficult, e.g. in Ref. |20| the authors provide an algorithm which does so but is 
exponential in the square root of the system size (see also L21J). 

To conclude this section, we study the possibility to use the data provided by the quantum experiments in 
order to construct a bound for the error on the estimated partition function. Our motivation is that, possibly, 
the a posteriori error analysis might be finer than the error bounds provided by plain error propagation. To 
simplify the discussion, we will again restrict ourselves to one-step protocols. Extension to the general case 
is straightforward. Let us expand the quantity A{i/3) as 

2N 

A{i/3) = - aj) + i ^w^^\i/3 - aj)) {^A{aj) + i ^A{aj)), (29) 
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Figure 2: Example of reconstructed quantities for a 10 x 10 classical Ising model with uniform ferromagnetic 
couplings (J = 1) and zero magnetic field (left) and for a 8 x 8 classical Ising model with non uniform 
couplings (J = ±1 with equal probability) and uniform magnetic field /i = 1 (right). The reconstructed 
quantities are the negative of free energy per spin (red), energy (green) and specific heat (blue) as a function 
of temperature and normalized by the number of spins. The plots show the average value of the quantities 
mentioned above, which is identical to the true value up to numerical machine precision, with error bars 
representing the a-priori standard deviation. For the simulation, we supposed to have experimental data 
with standard deviation equal to 10~^. 



and focus on, say, 

2N 

ARR{i^) = ^w^^\i^ - aj) ^A{aj). (30) 

The three other bits of A{i/3) are treated likewise. As was shown in the previous section, each quantity 
^A{aj) is obtained by measuring the polarisation of a qubit in a precise direction. Such a measurement 
process can be viewed as drawing a random variable whose outcomes are {+1,-1}, and whose mean value 
is the polarisation we are interested in. Let M denote the number of Bernoulli trials involved in determining 
each probability amplitude, and let us denote Xj{k) the outcome of the fc-th trial used in the determination 
of ^A{aj). Our estimate for ARR{i/3) is 

M 2N 
k=l j=0 

For fixed j, the random variables Xj (k) have the same distribution for all fc, characterised by Prob[Xj (k) = 
+1] = pj. If we assume there is no (uncontrolled) systematic error in the quantum experiments, then the 
true value of ARR{i/3) is of course given by 

2N 

ARRm = - aj){2pj - 1). (32) 

j=o 

Let us denote the second moment and third moment of a Bernouilli random variable X respectively as 
E2(X) = ((X-(X))2), E3(X) = (|X-(X)p). The following central limit theorem holds for the statistics 
of errors: 
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Theorem 4.1 (Central limit) Let E2(Xj) and E'^{Xj) denote (appropriate) estimates for E2(Xj) and E'^{Xj) 
respectively, obtained from M Bernoulli trials. Then, 

yS > 0, Prob[\ARR{i/3) - ARR{i/3)\ < S] > 
sup,>o{l - 2 6-^)^'^+'^ [1 - 2T,{-XM{e)S) - 1.12 Z)M(e)], (33) 



where 



mV2 

AM(e) = , ^ (34) 



DM{e) = 



1 j:^\^w(^\i(3-ajnUXj) + e) 



and where denotes the cumulative distribution of a zero-mean, unit-variance Gaussian probability dis- 
tribution. 

Proof: The proof of this result builds on the Berry-Esseen theorem ll22ll and is given in Appendix |Bj □ 



5 Link Invariants 

There exists several well-established connections between knot theory and statistical mechanics [23]. One 
of them is the following. For every knot it is possible to construct a graph such that the partition function of 
a Potts model defined on that graph is a link invariant for certain (imaginary) temperatures. This invariant 
turns out to be the Jones polynomial evaluated at specific values, modulo a known calculable factor. As 
the quantum algorithm for computing partition functions described in Section [2] is efficient for imaginary 
temperatures, it follows that it may also be used to distinguish among different link, when the associated 
statistical model only involves nearest neighbour interactions. In this section we outline the method to 
compute the statistical-mechanics knot invariant for any given link. We also compute these invariants for 
some primary linka with few crossings for which the Potts model involves only a few sites and is within 
reach of current technology. 

Let us start with a brief reminder on a recipe to construct statistical mechanical invariants, given a single 
component knot or a multicomponent link. We consider the planar projection of a given knot and shade the 
regions of the diagram in an alternating way such that there are no adjacent shaded regions (there are two 
ways to do this for any knot). We associate a lattice with vertices V and signed edges £, A = {E^ V) to the 
diagram in the following way. Every shaded region of the diagram will be a vertex of A and every crossing 
of the diagram that separates two shaded regions will be an edge linking the two vertices associated with 
those regions. The sign for the coupling of the edge is determined by the convention in Fig [3] For every 
edge z G f we associate a weight (a, a'), where a, a' are g-valued spins located at the vertices joined 
by the edge. Let us define a partition function given a set of weights Wi on L, 

{a} ieE 

where the sum is over all possible configurations of the spins on the vertices. 
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Figure 3: Convention for determining the sign of the edge coupHng assigned to each crossing (left). Example 
of the lattice obtained following the procedure outlined in the text for one of the possible shadings (right). 
Plain line represent, say, postive couplings, while dashed lines represent negative couplings. 
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Figure 4: Planar diagrams for some examples of primary knots and links (left) and the associated partition 
functions (right). Of the two possible graphs for each knot (one for each choice of shading) we have chosen 
the less trivial one. All the knots lead to statistical mechanics models with nearest neighbor interactions 
except for the Borromean ring, 62- 

Zl is invariant under ambient isotopy provided the weights Wi satisfy certain conditions, the derivation 
of which is discussed in |[23l . It has been proven that the choice Wf" = exp (±/35cr,cr^) where a = 1, . . . , g 
is compatible with these conditions if 



/3 cosh 



q-2 



(36) 



holds. In particular, the Potts partition function for g = 1, 2,and 3 at temperatures (3 = z27r/3, z7r/2, and 
z7r/3 respectively is a knot invariant. 

We have determined the lattices L for six examples of knots and links (see figure |5]) and computed Zl 
for a Potts model defined on L with g = 1, 2 and 3 for the values of /3 where the partition function is a knot 
invariant (see table[T]). The invariant corresponding to the value g = 1, 2 are actually trivial. The case where 
g = 3 is more interesting. A classical algorithm to compute this invariant exists which works in a time 
that scales polynomially with the number of crossings [24]. In turn, using a generalisation to three-level 
systems of the scheme presented in Sectionj2] allows to estimate the quantum invariant Zl in constant time 
with an additive error that scales like 1 / \/^where R is the number of repetitions of the experiment, now 
independent of the number of crossings. 
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Table 1: Knot invariants computed from the Potts model partition functions defined on the lattices in figure 
[5] The temperatures at which the partition functions have been evaluated are given in the text. 

6 Computational Power of Classical Models 

The analysis presented in Section|4]demonstrates how one can sample from a family of quantum circuits with 
fixed topology in d dimensions to construct a partition function on a classical spin system with fixed topology 
in d + 1 dimensions. One could ask whether the reverse can be done, i.e. given a classical partition function 
can one then reconstruct the outcomes of a related quantum circuit for a family of coupling parameters? 
Even more, is is possible that given the ability to compute the partition function of a suitably large classical 
system and for a suitable set of temperatures, one can reconstruct the outcome of measurements on arbitrary 
quantum computations of polynomial length in some fixed register input size? This has been partially 
answered in Ref. [2J where the authors show that the problem of computing the partition function of several 
classical spin models including the planar Ising model with magnetic fields all with complex couplings is 
BQP-complete. Such classical models do arise for some problems, e.g. the use of the Potts model with 
complex couplings to compute link invariants as discussed in Sec. |5| but it would be useful to have a 
result for physical, i.e. real, couplings. We do so in this section and also describe some applications: one 
for investigating quantum phase transitions given the ability to compute classical partition functions, and 
another for computing partition functions given the ability to prepare and measure corner magnetisation on 
physically prepared classical thermal states. 

6.1 Estimating quantum computations from Ising model partition functions 

We show the following: 

Theorem 6.1 Estimation of the partition function Z(/3) of a two dimensional ferromagnetic, consistent 
Ising model at inverse temperature f3 on a square lattice of size n x m with m — 0{poly{n)) with non 
uniform couplings and magnetic fields with additive error 6{n, m, /3) < exp(nm(49/3 — 190)/2) is BQP- 
hard, i.e. it is at least as hard as simulating an arbitrary polynomial time quantum algorithm on n qubits. 
By ferromagnetic we mean the couplings Jij in Eq. ^are all positive and by consistent the magnetic fields 
hi are all non-negative or all non-positive. To simulate a quantum algorithm means to do the following: For 
a unitary W built from a quantum circuit composed of 0{poly{n)) one and two qubit gates on a length n 
register provide an estimate of a complex scattering matrix element satisfying 

K+x |H/|+, ) |H/|+, 0{poly{n)) 

with a probability that is exponentially in n close to 1. 
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Proof: 

The proof follows in several stages. First we write an arbitrary polynomial sized quantum circuit in 
a convenient spatially translationally invariant form. Then we show that the scattering matrix element is 
equivalent to a complex temperature classical Ising model on a square lattice. Finally, we show that sampling 
the partition function over many real temperatures of a ferromagnetic Ising model, one can reconstruct the 
scattering matrix element. 

There are many possible equivalent quantum circuits which construct a given unitary. We pick a quantum 
circuit with a coupling graph given by a one dimensional chain of qubits with open boundaries. In order 
to perform a universal gate set, one needs a quantum circuit with gates either inhomogeneous in space or 
time or both. We pick circuits which are homogenous in space only as they are simple to parameterize 
and it is pedagogically satisfying that each step in the quantum algorithm can be thought of as a Wick 
rotated transfer matrix generated by a spatially homogenous quantum Hamiltonian. Several models exist 
for universal quantum computation which use ID architectures with global interactions ||T3l[T5|. We pick a 
convenient one due to Raussendorf |[T4l which involves encoding quantum information in a ID redundified 
data register, i.e. the data register is redundified in a second register which is spatially mirrored with respect 
to the first. This method has the advantage that all gates acting on the system are translationally invariant and 
the initial state is translationally invariant, e.g. |+f ^) . The only requirements are uniform Ising interactions 
between nearest neighbours and global single qubit gates. Addressability is afforded by temporal addressing 
via judiciously chosen homogenous local operations. Readout can be done again using global operations 
with the assistance of interspersed ancillary qubits or instead by using ancillary levels of each qubit ll25]| . 
The overall overhead incurred using global operations in this mirror encoded state is linear in n f\M. 

Consider a quantum register of an even number n of logical qubits, encoded by a chain of 2n qubits. 
The encoding has a mirror structure, i.e. the wave function of the system is at all times of the form 
|^)2n...n+i- The first ingredient in our proof of the BQP-hardness of the Ising model is the following lemma: 



Lemma 6.2 Let 



2n 



a. 



tot 



2n-l 



2n 



a^x.y, z, CPtot = W Hadtot = W Had^, (37) 



denote a set of translationally invariant ( global) operations, where Had = e 2^ ^ ^ ^ ^ denotes a single 
qubit Hadamard gate and CP — e^^'^^^^^^l the controlled phase gate. The subset 

e = {CPtot, <t(^/8), Hadtot} (38) 
is universal for quantum computation. 



Proof: This is proved in Appendix |Dj □ 
This lemma implies that for any 6 > 0, there exists a sequence of operators {Ct G©:t = 0...m — 1}, 



such that 

m—l 



i(+f"ii^i+r) - (+f "I n ^t\+fn\ < (39) 

where m = O (poly (log ^, ^)). Let (Jtot label classical configurations for the 2n-qubit chain (element of the 
computational basis). The action of (j^o^(7r/4) and CPtot (up to a global phase) can be expressed as 

<ot(^/8)ktot) = e'f-e^"^=^-^\a,,,), CPtotktot) = e^f ( -'^-fe+i) (40) 
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while the matrix elements of a Hadamard gate (up to a global phase) read 

(a|Had|aO = ^e^f (^+^Oe^f (41) 

These expressions will help us to express the quantum amplitude (+f Ht^o^ Ising 
partition function. It is convenient to introduce the following class of operators: 



tot 5 



where the exponents 60(5) and ei (5) take values in {0, 1}. Up to constant factors, it is clear that the operators 
^tot(^/8). Hadtot and CPtot can each be expressed either as a single T-type operator or as a product of at 
most 2 T operators. Consequently, we can write 

m-l M-1 

(+f "I n = "I n 'T'.i+f (42) 

t=0 s=0 

where M > 1. If M = 1 then the overlap is: (+f ^nj -Qm-i /:^|+|2n^ _ 2-''Zid{%) where Zi^^ is the 
partition function for a classical Ising model in ID with magnetic fields. Since one dimensional Ising models 
are exactly solvable for any temperature, including complex temperatures, then so is the overlap. Non exact 
estimations of scattering matrix element occur for M > 1. Since each layer operator Ct can be expressed 
as a product of at most two such operators Ts, we see that M is polynomial in n (since we assume that W 
is a polynomial depth quantum circuit). This last form of the quantum scattering amplitude, together with 
the identities ( |40l[4T] ) allow to express the quantum scattering amplitude as the partition function of an Ising 
model at imaginary temperature. Up to a global irrelevant phase, we have 

m—l 

(+f 1 n - ^ E^--"^^^ (43) 

t=0 {a} 

where H{a) denotes the Hamiltonian of the form ([T]), defined on a square (2n) x M lattice. Simple in- 
spection shows that all couplings (resp. fields) appearing in this Hamiltonian are positive integers, whose 
magnitude do not exceed 4 (resp. 17). 

Let us now assume we are provided with the following resource: 

IsingEstimator: Given an inverse temperature, /3, and an inhomogeneous Ising Hamiltonian, defined 
on a two-dimensional square lattice of size n^; x n^, a device provides an estimate Z{/3) for the partition 
function, Z{/3), that satisfies 

Prob[|Z(/3) - Z(/3)| < 6 5(n,,n^,/3)] > ^, (44) 

in a time that is polynomial in n^; , , /3 , 1 /e. 

Our goal now is to study how we could design the function S so that this resource allows for an efficient 
estimation of scattering amplitudes of quantum circuits. Since all magnetic fields and couplings appearing 



in the definition of the classical Hamiltonian associated with a quantum circuit are integers, the r.h.s of ( [43] ) 
can certainly be written as 



2n(M+2) ^ 2^(^+2) 

{a} k=-M 



15 



for some coefficients Ck- The value of the integer is at most max^- H{a). The r.h.s. of the last equation 
can equivalently be written as 



2n(M+2) 2^(^+2) ^ 

{^} 



where P is a degree- polynomial. For all /3 > 0, our resource allows to compute an estimate 

V{e-^) = e-^^'Z{p) for V{e-^) that obeys \V{e-^) - V{e-^)\ < ee'^^' 5{2n,M, p). Using a La- 
grange polynomial interpolation based on K points {{e~^^ ^V{e~^^))^j = . . .K — 1} (K > 2M' + 1), 
we re-construct the polynomial V as 



K 



z - e-^^ 



This reconstructed polynomial is in turn used to estimate our quantum amplitude as Ht^o^ ^^)est = 

^2n(M+2)^ 7^(6^^/-*-^). The error over this estimate can be bounded as 

m—l m—1 ^ K—1 

t=0 t=0 j=0 

It would be desirable to pick the integer K and the temperatures /3j in such a way that the r.h.s. of this last 
inequality is minimised. Presumably, calculus of variations might make this task doable. We have proceeded 
in a simpler way and made the choice 

e-^j = j/KJ = 0...K -1. 

Then, 

A closed form for the denominator on the r.h.s. of this expression can be easily worked out: 

J-l K-l 

in(^'-^)i=n(^'-^)x n {k-j)=j\iK-j-i)\ 

k^j k=0 k=j+l 

For the numerator, we observe that 

16^ 



TT |i^e-/i'5_^| ^ ^ TT |e-/i6_;t/i^| = exp [ V In J(cos - - ^' 

ill I |i^e-/i6_^-| 11 I / I |KeWi6_j-| ^LZ^ 16 



2 



sm 



k^j ' ' k=0 ' ' /c=0 

The argument of the last exponential is: 

K-l 

2 



Plugging these results in our bound for the error on the quantum amplitude, we find that, in the limit of large 



m—l m—l K—1 



|/,(8)2n| TT r I , (8)2n\ / , (8)2n| TT ^ i , (8)2n\ i ^ ^ (5(2n, M, (j/K) 6 

I 11 ^tnx 1 K^x I 11 /esti ^ 2ri(M+2) 1^6^16 _ ]U\iK - j - IV 

t=0 t=0 j=0 ' -'l-'-V J J' 
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Considering the case where K — 2M' + 1, having an error 

5(2n, M, P) < sin ^ e^P+^-^^^^^'2''^^+^')T((2M' + l)e-^ + 1)T((2M' + 1)(1 - e-^))(2M' + l)-^^' 
16 

is therefore sufficient for efficient reconstruction of quantum amplitudes. Note that the maximum energies 
from vertical and horizontal bonds in the lattice is 4(2n(M — 1) + (2n — 1)M) and the maximum local 
field energy is 2n{M - 2)8 + 2n8 + 2nM + 4M((2n - 2)2 + 2). Then we have the bound: M' < 
50nM — 12M — 24n. To work out how large 5(2n, M, /3) is compared to the partition function, we can 
compute the the needed accuracy for a function of the error 5(2n, M, /3') (we use a scaled temperature 
= /3/ In 2 to simplify the expression) 

ln^(2n,M,/30 ^ MMn2 , MMn 2(0.7387 + 2 log2(2M^ + 1) - 2 log^ MQ 



(3^2nM 2nM (3^nM 

For large system sizes, 

/error(2n> 1,M> > -25 In 2 + 50 In 2(2.7387)//3'. 
Finally we get a bound for the permissible additive error in the estimation: 

5(2n, M, 13) < exp(nM(49/3 - 190)) (45) 

Writing Ux = 2n, riy = M, since estimating Z{/3) for a polynomial number (linear in Uy) of temper- 
atures with additive error S{nx^ny^f3) on each provides the requisite estimate of the quantum scattering 
matrix element on a poly(n) sized quantum circuit, the complexity of the estimate of Z{/3) for an arbitrary 



temperature is BQP-hard. This completes the proof of Theorem[6[TJ □ 

We have found how much relative error we can tolerate in an estimation of a classical partition function 
and still accurately estimate quantum scattering amplitudes. How does this compare to known accuracy 
of classical algorithms which provide estimates of these partition functions? In Ref. |[26| Jerrum and Sin- 
clair construct a fully polynomial randomized approximation scheme (FPRAS) for computing the partition 
function of an arbitrary classical ferromagnetic Ising model that is consistent. Specifically they provide a 
classical algorithm that computes an estimate Z{/3) of the partition function Z{(3) = J2{cr} e^^^^^^J^^ for 
the ferromagnetic Hamiltonian H{{a}) on N spins, with a multiplicative error e and success probability 

3 



Prob 



\Z{^)-ZiP)\<eZiP) 



> 
- 4 



in a run time polynomial in N^l/e. This probability of success can be boosted to 1 — 5 in a num- 



ber log(l/5) of repetitions ll26l . Since the classical Hamiltonian in Eq. 43 is ferromagnetic, then when 
6{2n, M, /3) > Z(/3),IsingEstimator is no more powerful than FPRAS. In other words, if the tolerable 
error of IsingEstimator could be equal to or greater than Z{f3) for the relevant temperatures needed to 
reconstruct the scattering matrix element, then BQP-hard problems can be computed in polynomial time via 
FPRAS. This is not expected to be the case so we almost certainly have the requirement that the inequality 



in Eq. 6.1 is 5(2n, M, /3) < Z{^) over some significant range of temperatures and that it is smaller by 
an exponential in the problem size M^ Note it is known that the problem of exactly computing the par- 
tition function for even a ferromagnetic classical Ising model is #P-complete |26|. This complexity class 
is the same as that for counting the number of satisfying assignments of a Boolean function and counting 
optimal Traveling Salesman tours. Approximating the partition function with multiplicative error for an 
anti-ferromagnetic Ising model on a square lattice is NP-hard and for the ferromagnetic model but with 
general fields is approximation preserving reducible to the complexity class #BIS ETll . The latter is as hard 
as computing the number of independent sets, (an independent set is a set of vertices that does not contain 
both endpoints of any edge), in a bipartite graph which is thought to be of intermediate complexity between 
#P and FPRAS. 
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6.2 Ising models to compute quantum ground state overlaps 



We now consider an application of the mapping between classical partition functions and quantum scattering 
matrix applitudes: measuring ground state wavefunction overlaps of quantum Hamiltonians. It has been 
argued in Ref. [28] that wave function overlaps, termed fidelity overlap, can be a good witness to quantum 
phase transitions when the ground states straddle a phase transition point. In an ideal laboratory, this problem 
could be split in two: prepare two quantum registers in the desired states and measure the overlap using, 
for example, the protocols in Sec. [3] A possibility for the preparation step is to initialise the quantum 
system in the ground state |^^o) of some simple hamiltonian and to evolve this Hamiltonian to the target 
Hamiltonian . A fundamental result of quantum mechanics, known as the adiabatic theorem, is that if 
the Hamiltonian is modified slowly enough, the state obtained at the end of the evolution will be very close 
to the true ground state \G) ll29l . Crucially, the time of the evolution need only grow polynomially with the 
inverse of the minimum gap of the system, 7. 

The purpose of this section is to exhibit situations for which the adiabatic evolution need not be actually 
implemented. We are going to show that, in a precise sense, "time can be replaced with space". Roughly 
speaking, we are going to show that, instead of performing measurements on a quantum system of a given 
size, say "size" that has been evolved for a time "time", we can equivalently measure partition functions of 
classical Ising models prepared on a system of size 0(size x time). 

To make things precise, we will focus on the quantum transverse Ising model, described by the Hamil- 
toniarB 

where A denotes some d-dimensional lattice, and E{A) denotes the set of edges of A. We are going to view 
this Hamiltonian as a particular member of a family of time-dependent operators labelled by some time 
index, t. This family is 

H{t)^Ho + Hi{t), te[0:T] (47) 

where 

ieA {iJ)eE{A) ieA 

Without loss of generality, we will assume that /i^ > 0. In that case, |$o) = |+?'^') is of course the (unique) 
ground state of ^o- Evidently, H{T) = H^. The starting point of our construction is a discretisation of an 
adiabatic evolution 



Theorem 6.3 Let T satisfy the inequality 



\3 



T>r.(«,.)^i^?i™i™S:, (49) 

where 7 = i^ii^te[0:T]8^P H{t), where gap H(t) denotes the difference between the two lowest eigenvalues 
of H{t). Let L denote a positive integer, and let us define the discretisation step as 

T = T/L. (50) 

The quantity by which the state Ul-iUl-2 • ■ ■ ^^o|+?'^') deviates from the true ground state \G) of is 
at most 



/2(|/i| . lAI + iJiil . |£;(A)|) 
A = 5 + rW ^' ' ' ' I III I ^ + KL{\h\ . |A| + I J|| I . \E{A)\) • \h^\ . |A|T^ (51) 



^We use the hat on the operator to emphasise that this is a quantum Hamihonian. 
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where K is some constant. Each unitary Uk is defined as 

jj^^^-irHo^-irH,{kr)^ (52) 

This theorem, whose proof is given in Appendix |F| will help us to study fidelity overlaps, 

f = {G\G), 

where \G) is the ground state of and \G) is the ground state of some other Hamiltonian H^. For T{T') 
and L{V) large enough to build an approximation to / can be replaced in good approximation 

with 

/ W^WI . . . Wl_, Ul-iUl-2 ...Uo (53) 

For the transverse Ising model exemplified here the fidelity estimate which gives witness to a quantum phase 
transition is for the the case \G) being the ground state of with couplings /i = and \ G) being the ground 

state of with couplings = 0, J' = J, and h'_^ = h±-\- 6h±. Near the critical point, h± = J, there is a 
strong dip in the fidelity especially pronounced for 6±h/ J ^ 0.2 ll30ll . 

We are going to use a classical argument of quantum field theory |[3T1l , in a simple form adapted to our 
purposes, and show that the overlap ( [53] ) can be expressed as a partition function for a d + 1 -dimensional 
many body system at finite (complex) temperature, described by a suitable classical Ising Hamiltonian. The 
operator e~^^^^ can be expressed as the transfer matrix of a classical system, using the identity 1321 



r(/3) = ^e^^^V)(^'l + (54) 

Since on the other hand, 

e-'^^^^"" = cos(t/i^)(1 -ztan(T/i^)a^), 

it would be natural to make the identification e~^^ = —i tanr/i^, giving /3 = — ^ In tan (rh^), in order 
to relate ^rir2 to a classical model. Rather we are going to express the single-site unitary operator e"^'^^^^^ 
in terms of two operators T. For 6 > 0, let us define /3±(e) through 

e-^^±i^) =z^i(l±e). (55) 

One checks that 

r(/3+(6)) r(/3_(€)) = (2 - e2)g(/?+(e)+/3_(e))[i _ i^^a% (56) 

This choice of using two transfer matrices is not strictly necessary but it guarantees that the amount by which 
and /3_ need deviate from the imaginary axis is small which makes the connection to the traditional 

classical to quantum mappings lISTIl more transparent. 
So, for 

= tan(T/i^), (57) 
we see that the operator e~^^^^ can be expressed as a product of two classical Ising transfer matrices: 

|A| 



1-62 



e4 + 4 



This latter identity allows to express each operator Uk in terms of classical Ising transfer matrices. 
Introducing closure relations and bearing in mind that the operator Hi (t) is diagonal in computational basis, 
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Figure 5: Representation of the + 1 dimensional classical Ising spin lattice with couplings that encode 
information of the wavefunction overlap on a d dimensional quantum spin lattice A. Here the overlap is 
^^^1^^^ which is an approximation to the fidelity / = where \G) is the ground state of a Hamilto- 

nian H and \G) is the ground state of H^. The sequence Ul-i - - - Uq provides for adiabatic evolution, in 
small time steps r, of a time dependent Hamiltonian H{t) from the product state |+x '^') to 1^^"*") (which 
is an approximation to |G)), and similarly for the sequence Wl'-i • • • Wq, in steps which builds an ap- 
proximation l^^"*") of |G) from |+x '^'). Note that the number of gates L and to reach target ground states 
could differ as will the couplings generically. Each gate is a composition of diagonal gates with dimension- 
less coupling Hi{t) and two non diagonal gates with dimensionless couplings The temporal evolution 
of quantum gates can be represent on a classical spin lattice of one extra dimension with bond couplings as 
indicated on the left. For the Hamiltonian in Eq. [47]a>:(A:t) means dimensionless row couplings f3k between 
nearest neighbour spins, local fields of strength /3kh/ J, and couplings ;5±between rows. The parameters for 
adiabatic evolution to 1^^"*") are indicated with primes. 

the matrix elements of each operator Uk can now be expressed as a sum over paths on three copies of the 
lattice A: 

a{k) a{k+l) a{k+2) 

Here a{k) denotes a classical spin configurations over one copy of A, and the interaction jC, defined over a 
lattice A X A X A, is 

jeA jeA {ij)eE(A) jeA 

This interaction looks like a classical spin interaction with alternating complex couplings in the 

"time' direction which transfers between different copies of the lattice A and complex coupling within the 
lattice A. We would like to be able to chose variable couplings along the "space" and "time" directions so 
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we define a new interaction (assuming J ^ 0) 



'^a(fc),a(fc+i),a(fc+2) = P- EjeA + 1) + /3+ EjGA + + 2) 

A similar Hamiltonian can be written to represent evolution by gates : 

>^'<T(fc),<T(fc+l),<T(fc+2) = /3+ EjeA crj{k)crj{k + 1) + /3'_ EjeA + + 2) 

+/3'(L' + L - 1 - fc)(E(„)e£;(A) cTiik + 3)aj{k + 3) + ^ E,eA'^i(^ + 3))- 

Now we can write a Hamiltonian on the enlarged lattice A = {1, . . . , 2(L + L') + 1} x A, 

L-l L'+L-l 

cr(2/c+l),cr(2fc+2),cr(2/c+3) 5 

which takes exactly the form of a classical {d + 1) -dimensional Ising Hamiltonian but with complex cou- 
plings. The associated partition function depends on the vector of couplings /3 = /3, 13'} 

Z(/3) = J^e-^(W) 
and is a sum over classical configurations defined over A. 



Substituting this expression in Eq.(|53]), we see that the fidelity overlap can be approximated by 

1 r /l-e2li|A|r /I -e'2iL'|A| 



(58) 



where, e' is a solution to 2e'/ (2 — e'^) = tan(r'/iY) appropriate for the quantum Hamiltonian H* and the 
vector of complex variables /3* = /3* , is 



PI 



+ ^ log [l ± v^3 - cos(2t/i^) csc(t/i^) /\/2 - cot(T/i^))) j 



= T I f + ^ log ( 1 ± ( ^3 - cos(2t'/iY) csc(t'/i'^) / ^2 - coi^r'h'^ 

: 



(59) 



iJT 

^/2 , 



which were obtained by solving Eqs. 



55 



57 



Note for rh^ <Cl, = zbif zb ^ - Oir'^h?^). which is the 
statement that the analytic continuation is performed to nearly purely imaginary couplings strengths along 
the "time" direction. 

As described in Sec. |4] we can write the partition function in a power series in exponentials of the 
coupling parameters. For simplicity we assume J = J' = 1, /i' = /i = in which case: 



mi 1712 Tns m4 



2(«= E E E E E E 

gi=-mi g2 = -m2 g3 = -m3 g4 = -m4 g5 = -m5 gQ = -mQ 



r p/^9i p/3'g2 p/3+g3 p/3-g4 ^/^Lgs ^/^Lge 

^S'l, 5^2, 5^3, 5^4, 5^5 ,5^6^ e e e e ^ e , 



(60) 
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where: 



mi=m2 = L|A|, m3 = m4 = L'|A|, ms = L(L - 1)(|A| - l)/2, me = - 1)(|A| - l)/2. 

In Appendix |G] it is shown how the coefficients Cg^^g^^g^^g^^g^^g^ can be obtained by sampUng the partition 
functionfor 0(poiy(L^|A|, L'^|A|)) number of rea/couphng strengths (3 which then gives an estimate Z{I3). 
In order to obtain an estimate / of the fidelity overlap, we then need to perform an analytic continuation: 




1 -£2 nL|A| 



e'4 + 4 



L'\K\ 

z{n- (61) 



€4 + 4) 

Suppose we demand the error in the estimation of the fidelity to be e = 0{l/poly{\A\): 

\f-f\<e, 

and that the additive error in the estimation of the classical partition function satisfies 

Prob[\Z0)-Zm<e50)]>^. (62) 
Then it is shown in Appendix |G| that the following precision will suffice: 

60) < IQ{JT){J'T')L{L - l)L'{L' - i)|A|6e(/5++/5— 6-4)i^|Ale(/^;+/5--6-4)i>^|A| 

g(f-1.6)L(L-l)(|A|-l)g(| 



(f-1.6)L(L-l)(|A|-l)^(^-1.6)L^(L^-l)(|A|-l) ^^^^ 



To summarize, the required precision in the partition function estimation shows an exponential dependence 
on quadratic and cubic quantities in the system size. The origin of this dependence lies on the number of 
Fourier frequencies needed to reconstruct the partition function. This number is obtained by summing the 
amplitude of the bonds in the lattice (represented in Fig. |5]) associated with each coupling. The interactions 
corresponding to vertical (horizontal) bonds need a number of Fourier frequencies which is quadratic (cubic) 
in the system size. This different behaviour ultimately comes from the chosen adiabatic time dependence on 
the total Hamiltonian of the system (Eq. [47| ). 

We close this section by commenting that while fidelity overlaps could be estimated using the method 
of mapping to a generic quantum circuit presented in Sec. |6.1[ the method described here is much more 
efficient in resource scaling since the gates are applied directly using the transfer matrix formalism rather 
than mapping to a fixed library of quantum gate in an encoded circuit. Furthermore, the required accuracy 



of estimation of the partition function is exponentially better than the bound computed in that case (Eq. 45 ) 



6.3 Corner magnetisation and estimating partition functions 

The foregoing analysis illustrates the computational power of accurate evaluation of Ising partition func- 
tions. We can wonder what is the computational power of more modest tasks, such as estimating the mean 
values of specific observables. We have studied a simple instance of this problem. As it turns out, very 
simple tasks already have computional power. For instance, the ability to accurately estimate single site 
magnetisations on random Ising models lead to random approximation schemes for partition functions. 
This is the content of the following theorem. 

Theorem 6.4 Consider the Ising model on a two-dimensional square lattice A, described by the Hamilto- 
nian: 

H{a) = -J^^aiaj -h^^ai. (64) 

{ij) ieA 
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For any e, inverse temperature /3, and magnetic field strength h it is possible to provide an estimate Z{^, h) 
for the Ising partition function Z{/3^h) satisfying 

Prob[|Z(/3,/i) < eZ(/3,/i)] > 3/4, (65) 

in a time that scales at most polynomially with l3,e~^, \h\, and the size of the system if we are able to 
perform corner magnetisation measurements on specific non-homogeneous Ising systems with a relative 
precision that need not be lower than the inverse of some polynomial in \h\, and the size of the system. 

Proof: The proof is given in Appendix |E| □ 

This result might appear surprising since it applies to even to antiferromagnetic Ising models whereas, as 
discussed above, a multiplicative approximation of the partition function in that case is an NP-hard problem. 
However, corner measurement is a quantum process which assumes the thermal state of the classical Hamil- 
tonian has been prepared. Some recent work ||20l[34l provides quantum algorithms to simulate thermal states 
of classical spin models. However as mentioned in Sec. |4j generically these algorithms scale exponentially 
in the system size, and given the complexity of multiplicative approximations of antiferromagnetic partition 
functions we would not expect a drastic improvement in thermal state preparation by quantum algorithms in 
that case. Whether efficient quantum algorithms exist for preparing ferromagnetic thermal states is as far as 
we know an open problem but if so than corner magnetisation measurement could prove a useful diagnostic 
for such algorithms since classical FPRAS is available. Finally, we add that recently quantum algorithms 
for FPRAS were found which exhibit a quadratic speed up over the classical counterparts |35|. These al- 
gorithms are rather different in spirit from measuring corner magnetisation as instead of using mixed states 
they use a combination of Grover search and phase estimation to prepare pure states of many qubit systems 
which coherently encode probability distributions of various classical spin configurations. 



7 Conclusions 

In conclusion, we have presented schemes allowing for the measurement of partition functions and mean 
values of classical many-body systems, at complex temperatures. Although we have mainly focused on 
Ising Hamiltonians, these schemes can be generalised to other systems, such as the g'-state Potts model for 
instance. We have presented two applications of these schemes. 

First, we have studied the possibility to use it in order to compute real temperature partition functions. 
Although our findings yielded results as poor as previous attempts made by other authors, it is interesting 
to have found similar results using a different route, in particular one that involves reconstructing partition 
functions for all temperatures as opposed to a single temperature. We have also seen how experimental data 
allow to a posteriori sharpen error estimates, through a central-limit theorem. To the best of our knowledge, 
the problem of determining whether quantum mechanics can be used (or not) to efficiently compute partition 
functions of classical models, or even FPRAS thereof, is still open. As a second application, we have 
seen how some link invariants could be deduced from the ability to detect imaginary temperature partition 
functions, again using constant depth quantum circuits. 

These applications all rely on two kinds of schemes, one whose implementation could, in principle, 
only require a constant time, another involving a time evolution. All schemes translate naturally into global 
operations and measurements supplemented by edge addressability. This is natural for certain architectures 
such as cold trapped atoms in optical lattices |36|, or superconducting qubit arrays |37|. Furthermore, this 
kind of quantum processing can be made fault tolerant without demanding more addressability as shown in 
1381. 
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We have considered the dual of the first application mentioned, and studied the possibility to efficiently 
simulate a quantum computer, given the ability to estimate real temperature disordered Ising partition func- 
tions. We have found that quantum amplitude of a depth- quantum circuit, acting on n qubits, could be reli- 
ably estimated if suitably associated disordered Ising models could evaluated with a precision that essentially 
grows exponentially with D and n. The problem of simulating quantum circuits from statistical mechanical 
partition functions, estimated with a looser precision (polynomial, say) is, just as open its dual. One implica- 
tion is that given the power to compute classical partition functions in + 1 dimensions, in certain cases one 
can compute quantities relevant to quantum phase transitions in d dimensions. This argument involved view- 
ing the overlap of two ground states of a quantum Hamiltonian as the scattering matrix element for a quantum 
computation which can then be estimated by computing classical Ising model partition functions with real 
couplings. The method was illustrated for the particular case of the quantum transverse Ising model in one 
dimension and while that model , the mapping extends to a variety of other quantum spin Hamiltonians in 
a straightforward manner. For the sake of simplicity, we have restricted our analysis to some homogeneous 
quantum Ising models in d dimensions and that model already has a well known correspondence to the 
(i + 1 dimensional classical Ising model |[23]| . However, extensions of our ideas to disordered quantum spin 
Hamiltonians of the form, say, H = (K^f + + K^i) " E{i,j> (Jij^f'^j + ^Ij^^i"^] + ^ij'^i^j^ 
is straightforward]^ This mapping could provide new ways to perform quantum simulation, via either quan- 
tum or classical algorithms for estimating Ising model partition functions. Given some of the difficulties that 
beset fault tolerant implementations of quantum simulations ||39l|40^^| new approaches are certainly desirable. 

Finally, we have seen how the ability to prepare thermal states and perform single qubit measurements 
immediately implies the existence random approximation schemes. This observation naturally leads to 
wonder what is the quantum complexity of the preparation of classical thermal state. In view of recent 
inapproximability results 11411 . it would be very interesting to solve this question in the case of the anti- 
ferromagnetic Ising model for instance. 
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A Disordered Systems 

Preliminary: We found it convenient to use a slight variation of the detection schemes described in Section 
[2] and consider single qubit gates defined from two Hadamard gates and a single phase gate: 

G{9) = Had ( J J, ) Had . (66) 

For 0* = —i log tanh /3 J, this single qubit gate turns out to be equal to r(/3 J)/2 cosh(/3 J), where T(/3J) 
is the two-spin Ising transfer matrix introduced in Eq.([54|). 



^Again, one could use the Baker-Campbell to decompose the evolution operator associated with this Hamiltonian. Then, it 
would be enough to express cr^ in terms of and operators using an Euler angles decomposition. 
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In this appendix, we are interested in two-dimensional Ising models, of size n x m, with random bond 
interactions, which strengths take value in { — 1, +1}. The magnetic fields felt by each spin is also assumed 
to be random and takes value in { — 1,0,+1}. For fixed configuration of bonds and magnetic fields, the 
partition function can be evaluated for a specific range of complex temperature, with the help of a two- 
dimensional lattice of quantum particles on which instantaneous measurements are made, or through the 
time evolution of a one-dimensional quantum system. 

The one-step protocol doesn't pose any particular problem for disordered systems. From quantum am- 
plitudes of the form given by Eq.Q evaluated at specific angles, one can reconstruct the partition function 
through analytic continuation. Namely, 

= A(z/3) = J2 ^^"^'^(^/^ - ^jM(jf)^ (67) 

ji=0 ^ 



where w^^^"^ is defined by Eq.(26), and where A^i is polynomial in n and in m. 

The case of the time evolved scheme is slightly more complicated than in Section [2j Reproducing 
the reasoning presented in that section, one can find an appropriate sequence of controlled phase gates ^ 
and G— gates that provides relevant quantum amplitudes. The real partition functions are again obtained 
after Fourier transform and analytic continuation. It turns out that three parameters are enough for that. 
One, a, takes into account constant-time interactions and magnetic fields. The two others, 9^ and 9~, are 
respectively related to ferromagnetic and antiferromagnetic interactions between particles corresponding to 
consecutive time-slices. More precisely, one can see that the kind of partition functions we wish to consider 
can be written as 

where A^i , A^2^ , N2 are again polynomial in n and in m. Actually, Ni = 2nm — n represents a bound on the 
total number of "horizontal" bonds plus the number of sites, while A^2^ (resp. A^2^) represents the number 
of ferromagnetic (resp. antiferromagnetic) edges connecting spins at different time-slices ("vertical" bonds) 
(A^2^ + A^2^ = m(n — 1)). The coefficients c^^^+^- are essentially Fourier transforms of quantum amplitudes 
A{a, 9^,9~) detected at selected angles a, 9^,9~ E (0, 2ti]\ 

Nt N- + ^ + 

c - = i-lY^ V V V — A( (69) 

^1^2+^2 ^ ^ (2iVi + l)(iV2^ + l)(7V2" + l) 2iVi + l'7V2^ + l'iV2" + l 

32 ~^32 ~^ 



One can note how the particular form of the G— gate (which does not involve terms of the form e ^^) allows 



for the Fourier transform in 9^ to be restricted to positive frequencies. Again, plugging Eq.(69 ) into Eq.(68) 
allows to express the partition function as a function of the "experimental" data: 

2^(e^ + e-^)^2++^2" 

^^^^ ^ (2iViH-i)(iv^ + i)(A^,- + i)S .f ^ 

2l7TJl _2VKJ^ 



(e/3e"^^)-^i 5(2^1 )(e^e"^^)5(^)(tanh/3e ^2++i tanh^e ^2+1), (70) 
where S^^\q) = (1 - q^+^)/{l - q). 
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The restricted set of possible values for the couplings and magnetic fields implies that the partition 
function of the disordered Ising model we are considering can be written as 



N 

m = ike-^^ (71) 

k=-N 

where again N scales polynomially with the system size, and where each ^/^ is positive integer whose 
magnitude is at most equal to the number of possible configurations for the system, i.e. ^ 2^^, V/c, i.e. 
they can be represented exactly with nm bits. Thus, the estimation of each coefficient with nm bits of 
accuracy, i.e. with a variance ^2{ik) lower than one would allow for an exact reconstruction of the partition 
function for all temperature. Yet another Fourier transform allows to show that 

27V + 1 ^ ^ 27V + r ^ ^ 

j=0 



Combining this latter relation with Eq.([67|) for instance, it is possible to see that in order to get with nm 
bits of accuracy, one would need to estimate the quantum amplitudes themselves with 0{nm) bits of accu- 
racy. Unfortunately, we do not know how to do that efficiently. In our scheme, the quantum amplitudes are 
obtained from repeated Bernoulli trials. It therefore seems that 0(2^^) trials are then necessary. A similar 
conclusion is reached when the evolution protocol is used. 

We now give some more qualitative insight on the performance of the protocol by analizing a particular 
instance of the reconstruction (through the time evolving algorithm) of the coefficients (Eq. U^ for an 
8x8 Ising model with 50% positive/negative bonds and uniform magnetic field (set to 1). In Fig. [6]we plot 
the coefficients togheter with an upper bound on their standard deviation as a function of k. 
We remark that, for small temperatures, only coefficients with big k are important as it is evident from 
the series in Eq. ITT] As shown in the plot, in the "big fc" range, two facts are evident: the standard devia- 



tion goes to zero and the coefficients are exactly zero. The reason behind the behaviour of the standard 



deviation is found by algebrically expanding equation [68] and noticing that the coefficients (responable for 
the amplification the experimental errors) multiplying big powers of are small. On the other hand, the 
behaviour of the coefficients for big /c is a natural feature of the disorderd model consider here. More 
specifically, it simply reflects the impossibility for the ground state spin configuration to minimize each lo- 
cal term of the hamiltonian, namely, to satisfy each bond and align to the magnetic field everywhere. The 
low temperature properties of the model are then shifted around the fcs where the coefficients start to be 
non-zero. Unfortunately, in that regime, the error is not longer approacing zero, explaining why, in this case, 
the protocol does not give good low temperature performances. Conversely, for a uniform ising model, the 
coefficient would be nonzero for the biggest possible k explaining how, in the main text (see Fig. [2]) we 
could get good results in the low temperature limit for this uniform case. 



Again by simply inspecting Eq. [71] one can easily get convinced that the importance of the coefficients 
for smaller k becomes more important as the temperature increases. In this regime, the standard deviation 
is basically constant witnessing properties of the counting process needed to calculate the coefficients 



again obtained by expanding eql68|in powers of e^. As one can infer by the plot, the relative error is quite 



small for these intermediate values of k, justifying the better high temperatures performances of the protocol. 



At last, we would like to qualitatively analyze the possibility to use the protocol to estimate an upper 
bound of the ground state energy. This relies on restating the problem of finding the ground state energy 
as the problem of finding the maximum k for which ^^^^^ ^ 0. Following this statement, in order to find an 
estimate for the upper bound for the ground state energy, we want, roughly speaking, to look at the condition 
by which the standard deviation on the coefficients is not bigger than the coefficients themselves. In the 
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Figure 6: (color online) Plot of the coefficients (in blue, see Eq. |7T]) and an upper bound on their 
standard deviation cr^^ (in red) as a function of A: for an 8 x 8 Ising model with 50% positive/negative bonds 
and uniform magnetic field (set to 1). From this plot we can qualitatively justify the performances of the 
algorithm in the small and high temperature limits. The low temperature limit behaviour has to be found in 
the range of k where the coefficients start to be non-zero. This range does not correspond to the maximum 
possible value of k due to the fact that, in the present model, spin configurations cannot minimize each local 
hamiltonian. This does not allow to take advantage of the enhanced precision of the protocol for big k 
and it is the reason for the poor performances of the algorithm at small temperatures. As the temperature 
increases, the whole range of k starts to become important, so that we can focus on the intermediate values 
of fc, where the bigger coefficients As evident from the plot, in this regime the relative error is 

quite small explaining the good high temperatures performances of the protocol. The value of k where the 
standard deviation is equal to the relative coefficient sets the limit for a possible estimate of an upper 
bound on the ground state energy. 
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Figure 7: Logarithm of the relative errors in the reconstruction and analytic continuation of the square of the 
partition function of the classical ferromagnetic Ising model with open boundary conditions as a function of 
the inverse temperature. Each curve is the relative error for a different system size. For all system sizes the 
error shows a peak near the critical point. However, for systems larger than 7x7 the error grows quickly as 
the temperature approaches zero. 
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plot presented here, an upper bound on the ground state energy is then obtained by looking at the point 
where the two curves intersect. As one can see, the result for this instance is very good, but, generically 
speaking, the impossibility to rule out worst cases scenarios, does not allow us to give more quantitative 
results. 

In the same way that protocol 2 can be used to estimate the real temperature partition function via mea- 
surements of the quantum overlap {^l"^) and an analytic continuation, protocol 1 can be used to estimate 
the square of the real temperature partition function via measurements of the square of the overlap | ($|^^) p. 
The only difference is that as the function to reconstruct is squared, the frequencies of the modes in the 
Fourier series that we construct from experimental data is doubled. Hence, for protocol 1 more measure- 
ments are needed, double the amount needed in protocol 2. We have reconstructed the square of the partition 
function of a classical Ising model in 2D and performed the analytic continuation. For a study of the errors 
in the reconstruction see FigJTJ 



B Proof of theorem [13] 

Let us start with the following classical result |[22ll 



Theorem B.l (Berry-Esseen) Let Wi . . . Wl denote L independent random variables such that (Wk) = 
0, < (VF^) < oc, < oc, k G {1...L}. The cumulative distribution function Tw of 



(a? + ... + 4)1/2 
satisfies the inequality 



\\:fw-:fa\oo < Cbe Y^^Wu?) / {Y.^wl)f'\ m 



.2\^3/2 

vv 

k=l k=l 



where denotes the cumulative distribution of a zero-mean unit-variance Gaussian. The value of the 
constant Cbe is most 0.56 / l?2l/ . 

In the context of the estimation problem of Section|4j the r.h.s. of Eq.(|73]) reads 

Cbe Dm = Cbe2_^ — j . (74) 

The quantities E2{Xj) and Es(Xj) are of course defined as in Section]^ Their values cannot be known 
exactly. But, from experimental data, one can construct estimates for them, E2 {Xj ) and E3 {Xj ) respectively, 
valid with high confidence. With such estimates, we can bound Dm as 

Dm < DM{{e}) = ^ ^ | W^)(i/3-a,)P(E3(X,)+e,)/( ^ |5ftu;W(i/3-a,)|2(E2(^^^ e, > 0, 

J J 

at a confidence level YljFYoh[{E2{Xj) > ^(Xj) - e^} n {Es(Xj) < Es(Xj) + e^}]. ({e} is a short- 
hand notation for {ej}0. Let dA denote the random variable Ajiji(i/3) — ARji(i/3), and let J'dA denote its 
cumulative. 
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By definition, Probflrf^l < 5] = J^dA{^) — ^dA{—^)' As just any cumulative distribution, J^^A satisfies 
the scaling identity TdA{y) = J^x (iA(Ay), VA > 0. Thus, with 



Am = ^M/^|3?«;W(i/3-a,)|2E2(X,), 

we can use the Berry-Esseen theorem to lower bound the probability of the event {\dA\ < S}, conditioned 
on the event {Dm < DmH^})}- 

Froh[{\dA\ < 6}\{Dm < DM{{e})}] > 1 - 2J'*(-Am5) - 2Cbe DM{{e}). 

Of course, cumulative functions are non-decreasing. Thus, 

Fvoh[{\dA\ < 6} I {[Dm < DM{{e})] H [Am({6}) < Am]}] > 1 - 2J'*(-Am({6}) 6) - 2CBEDM{{e}) 
where 

AM({e}) = ^M/ |5ft^i;W(i/3 - a.WiUXj) + e,). (76) 

Finally, using Bayes' rule, we find that 

Pwh[\dA\ <S]^ Pioh[{\dA\ <S}\{[Dm< DM{{e})] n [AM({e}) < Am]}] x 

Prob[{DM < DMi{e})} n {AM({e}) < Am}] 
> nProb[{E2(X,) + > E2(X,) > E2(X,) - e,} n {E3(X,) < E^iXj) + e,}] x 



J 



[1 - 2J'*(-Am({6})5) - 2CBEDM{{e})] , V5, {e} > 0. (77) 
We show in Appendix |C] that 

Prob[{E2(X^-)+^j > ^2{Xj) > E2(X^-)-^j}n{E3(X^-) < E3(X^-)+^j}] > 1-2 exp(-e^-M/32), Ve^- > 0, 
So, in particular, for fixed S, 

FYoh[\dA\ <6]> supj]^ (1 - 2exp(-e^M/32)) [l - 2J'*(-AM({e})5) - 1.12 DmIW)], (78) 



from which Eg. ([33]) follows. 



C Confidence intervals for moments of a binomial 

Consider a binary random variable V characterised by some probability p for one of the values. Let p denote 
our estimate for p, found after M Bernoulli trials. Let E2(p) denote the variance of V. Clearly, 

\p-p\<eo^ E2{p) G [E^(p,6o),E+(p,6o)], (79) 

where E2"(p, eo) = min{E2(g) : \q — p\ < ^o} and E^(p, eo) = max{E2(g) ' \q — p\ < ^o}. Similarly, if 
we denote E3(p) = (|y p), we have that 

\p-p\<eo^ Es{p) G [E-(^,6o),E+(^,6o)], (80) 
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where Eg (p, eo) = min{E3(g) : \q - p\ < eo} and E^(p, eo) = max{E3(g) : \q - p\ < cq}. For 
bookkeeping, let us attribute names to the various events we need to consider: 

A = {\p-p\<eo}, B = {E2{p) e [E2-(p,eo),E+(p,eo)]}, C = {E^ip) G [E^ip,eo),Et{p,eo)]}. 

Clearly, A c B and A c C, so that ^ C -B n C, and 

Prob[S n C] > Prob[A] > 1 - 2exp(-€oM/4). 



The last inequality is a consequence of Hoeffding's II221I . Next we define our estimate for the variance and 
the third moment as 

2 _ E2 (P, ep) + ^tiP^^o) _ E3 (p, ep) + ^tiP^ gp) 
' ~ 2 2 

We also define 2ei = E^(p, eo) — E2"(p, eo) and 2e2 = ^s{p^ ^o) — E3 (j3, eo). For any e greater or equal to 
max{ei, 62}, we have B d B' and C C C^ where the events and C are defined as 

B' = E2(p) G - 6, 5^ + e], = Eab) G [r - e, r + e]. 

So, Prob[S' n > Prob[S n C] > 1 - 2 exp(-eoM/4). Picking 

dE2(2?)| dE3(j)) 
6 — maxj sup — |g, sup — \q} x eg = eg, 

0<g<l dp 0<q<l dp 

we find that 

Prob[{5^ - e < E2{p) < 5^ + e} H {E3 < r + e}] > 
Prob[{52 - e < E2(p) < 5^ + e} H {r - e < E3 < r + e}] > 1 - 2exp(-e M/32). (81) 

D Proof of Lemma \6.2\ 

We begin with the discrete gate set 

©0 = {Zfc(V4), Hadfc, A: = 1 . . . n } U {CNOT^^^+i, A: = 1 ... n - 1} (82) 

acting on an n qubit register that is universal for quantum computation ll43l . By the Solovay-Kitaev ll44l 
theorem an arbitrary polynomial sized quantum circuit can be efficiently approximate from this gate set with 
a polynomial overhead. To realize this using global operations in the mirror encoding of Raussendorf, one 
makes frequent use of the global shift operator Gtot = a^Q^{7v)afQ^{7v/2)CPtou with the property that G^^^^^ 
is a reflection of the state of the chain about its middle. An arbitrary Z rotation on logical qubit k can be 
physically implemented as f\M 



(83) 

Similarly, an X rotation on logical qubit k is 

X'f '(a) = e^f = Gl-^ aM G al,{n) G^ af„,(-7r/2) aUa/2) al,{^/2) 

Gro7' C(vr)Gc7f„t(7r)G'=c7f„,(-7r/2) c7f„t(-a/2X,(7r/2). (84) 
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Finally, an entangling gate between logical qubits k and A: + 1 can be implemented as 

Vg^,(a) = e^"("^^"^^+i+"^+-^"^+-i) = G^X'^^\a)G^\ (85) 

Since T4,/e+i(vr/4)Hadfc+i^/e(vr/2)Zfc+i(7r/2)Hadfc+i = CNOTfc,fc+i then the gate set 

©1 = {Z'^'\7v/4), HBd'^^\ A: = 1 . . . n } U {Vk,k+i{7v/4), A: = 1 ... n - 1} (86) 

is universal for quantum computation. Now the Hadamard gate can be related to X and Z rotations through 
the identity Had = a^(7r/2)a^(7r/2)cr^(7r/2). Also we note the following relations: [cr^(7r/8)]^^ = 
a%-7T/8),mday{±7T/4) = a'^{-7T/2)a^{^7T/4)a''{7T/2) and also ^^(±^2) = a^{±7T/2)Ha6a%±7T/2). 
Then from Eqs.( |83|84|85| ), we see that it is enough to be able to implement 

6 = {CPtot,af,t(V8),Hadtot} 
in order to achieve universal quantum computation. 

E Magnetisation and Approximation Schemes 

The first part of our construction closely follows a general argument presented in Ref. EH, and establishes a 
connection between partition function evaluations and the ability to draw samples from Boltzmann probabil- 
ity distributions. Some adaptations were made, though. We felt that indicating only these adaptations would 
have resulted in an awkward presentation. This is why, for the sake of clarity, we have chosen to reproduce 
this argument, with these adaptations included, in a concise but self-contained manner. In the second part of 
our construction, we show how measurements of magnetisation on specific non-homogeneous Ising models 
allow to draw from Boltzmann distributions. Consider thus the Ising model on a two-dimensional square 
lattice A, described by the Hamiltonian: 

H{a) = -J^aiaj -h^ai. (87) 

{ij) ieA 

For /i = 0, the model is solvable and Z(h = 0) is known exactly (see e.g. |[23l ). We wish to evaluate the 
partition function at a fixed temperature|j/3, Z{h), for h > 0, sa}j^ For that purpose, we express Z{h) as 

where = ho < hi < . . . < Hl = h. These values hk are chosen to be equally spaced, and we will denote 
the spacing h^ — h^-i by Sh. Each ratio = Z{hk) /Z{hk-i) can be expressed as 

Qk-Y. ^^^e/^^^l^l^(-) ^ e/^^^l^l^(-), (89) 

where M(cr) denotes the mean magnetisation of the system when the lattice is in configuration a, |A| 
denotes again the size of the lattice A, and where iJ^-i is a shorthand notation for the hamiltonian when the 
magnetic field is set to h^-i. 



"^Change of notations: Since we will work at constant temperature, we will from now drop /3 and simply write Z{h) instead of 
Z{p,h). 

'The case ft < is treated similarly. 
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In order to evaluate we will use a collection of estimators for the quantities Qk, each involving n 

sample configurations. These estimators are defined as 

: at'} ^.(a^K • - - , 4"^) = ^ E e^^^S.^i^i'), (90) 

where each sample a^^^ is drawn according to some probability distribution 7r^_^. Our estimator for Z{h) 
is 

L 



k=i 

Let Qk denote the mean value of Qk, i.e. 

Since all ^^^^^ are independent random variables, we find that the mean value of Z{h) is given by Z{h) = 
Ylk=i QkZ{ho). Now let us assume that 

\Z{h)-Z{h)\ < eZ{h), (91) 

and that 

\Z{h)-Z{h)\ <S Z{h), (92) 

with probability at least, 3/4 say[^ Then 

(1 - S){1 - e')Z{h) < Z{h) < (1 + 5){l + e')Z{h), 
with probability at least 3/4. Thus 

(1 - e)Z{h) < Z{h) < (1 + e)Z{h) (93) 

with probability at least 3/4 whenever e > 5 + + 5e^ 
Clearly, 

Plugging these inequalities into Hoeff ding's inequality [22], we find that 

Prob[|^fc -gk\<C]>l- 2e-2<'/^i^M|A|/3c5/^)2^ ^94^ 

Let us use this latter relation in order to construct an upper bound on \Z{h) — Z{h)\ valid with tunable 
probability. We will use the following Lemma: 

Lemma E.l 

\Z{h) - Z{h)\ < I 17(1 + ^) - 1| Z{h) (95) 

Li 

with probability at least (1 - 2e-2<Vsinh(|A|^5/z)2^L 

^This value is somewhat arbitrary. As explained in Ref.||26|, any level of confidence strictly above 1/2 can be efficiently boosted 
to arbitrarily close to 1. 
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Proof: We start with the following identity 

L L 



\z{h) - z{h)\ = I n - n ^(^o) = 1 11(1 + ^^^) - ii z{h) 
k=i k=i k=i 

Next, we have the inequality 

L L 

I + - ii < I + i^^i) - ^ ^- (96) 

k=l k=l 

Let us consider two cases: (i) 11^=1(1 + ^k) — 1^0, (ii) 11^=1(1 + ^k) — 1 < 0. The inequality is trivial 
in case (i). In case (ii), we need to prove that 



l_Y[{l + Xk)<l[{l + \xk\)-h 



k=l k=l 

or 2 < 11^=1(1 + \^k\) + 11^=1(1 + ^k)- The r.h.s. of this last inequality can certainly be written as 



L ^ 



where each coefficient x^^...^^ is non-negative. It is also clear that each quantity . . . . . . 

is non-negative. Inequality ( |96| ) implies that 

Qk - Qk\ 



\Z{h) - Z{h)\ < I 11(1 + IMK^m) _ i| z{h). 



k=i 

Ther.h.softhisrelationislowerthanI nLi(l+i)-l| ^(^) with probability at least (l-2e-2<Vsinh(|A|/35/i)2>jL 
(Hoeffding's inequality). □ 

We will pick the spacing between two consecutive magnetisations to be Sh = ^j^, where 77 is some 
positive constant we are free to choose at our convenience. Sh fixes the value of L to 

L = (/i-/io)/3|A|/77. (97) 

With a given choice for 6h, we have that Qk > and 

\Z{h) - Z{h)\ < 1(1 + e^C)"^ - 1| Z(h), (98) 

with probability at least (1 - 2e~^^^^/^^^^(l^l^'^^)^)^. How should we pick ( in order to ensure that the l.h.s. 
of (98 ) is smaller than 6Z(h) for some fixed 6 ? Since (1 + e^C)^ < e^''''^ it is enougli tliat 

^ - Lev ' 

We also wisli to know how, for fixed values of (, L, r], we should choose n in order to guarantee a level of 
confidence at least equal to 3/4. Direct substitution shows that the condition 

(1 - 2e-2"^'/''°''('')')^ > 3/4 
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is satisfied if 

sinhTi^e^^L^ ^ 
"^- 2(ln(l + ^))2 ^" 



(99) 

To summarise, for L satisfying ( [97] ) and n satisfying ( [99] ), the partition function estimator satisfies 

Prob[|Z(/i) - Z{h)\ < SZ{h)] > 3/4. (100) 
Next we wish to estabhsh a condition that guarantees that InequaUty ( [9T] > holds. We start by observing that 

\z{h) - z{h)\ < 1 17(1 + l^hJi^) - 1| z{h). 

Let ^ 

ATTfc-i = max|7rfc_i(5) -4-i(5')| = - |7rfc_i(cj) - 4_i(cr)| 

cr 

denote the total variatioij^ between the probability distributions Hk-i and 7r^_^. Let us also denote Att* = 
max/e A7rj[;^_i. We see that 

Thus 

- < [(1 + e^^A^*)^ - l]Z(/i) < (e^^'"^^* - (101) 



So it is enough that 
On another hand, Ank-i satisfies the inequality 



A^* < — — ln(l + eO. 

1j 



Att^;-! < - max 1 — . 

2 ^7 7r/e-l(cr) 

Using Bayes' theorem, to express tt^-i in terms of marginal and conditional probability distributions, 

7Tk-l{(Ji . . .cr|A|) = 4-l(^l)4^-l(^2|cri) . . . ^'^f (a|A| |cri . . .cr|A|_i), (102) 

the r.h.s of the latter inequality can be written as 

- max --r^ 1 . 

If we use the finesse 

J = max (103) 

to quantify the accuracy with which the distributions {7r^_^} approach the distributions {nk-i}, we see that 
Att* < ^1(1 + f)l^l - 1| < ^(e^l^l - 1). So \Z{h) - Z{h)\ < e'Z{h) whenever the finesse satisfies 

1 2p~^^ 

f<^ln[l + ^ln(l + e')]- (104) 
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Figure 8: A. Typical corner on which magnetisations need to be measured in order to get an approximation 
for the partition function of the Ising model defined on a square lattice. B. Labelling of spins of the original 
lattice. C. Lattice obtained after the spin on one corner has been fixed. D. Cartoon for a possible choice to 
run over all spins of the original lattice 



We now turn to the second part of our construction and explain how it is possible to get samples for the 
estimators Qk from corner single site magnetisation estimates, as indicated on Figjsj-A. Assume that the | A| 
particles of the lattice are numbered as indicated on Figj8]-B. For fixed external field hk-i. It is clear that the 
magnetisation on the corner 'V is given by 

From an estimate m'^_^(l), we construct 7r^^l^-j^((Ji) as 

4i(t) = ii%^. = (105) 



It is certainly possible to draw exactly according to this distribution vr^^^p it is a known two-outcome prob- 
ability distribution. Let us imagine we do it and obtain an outcome xi. Then we consider another Ising 
system, identical to the original apart from the fact that the spin labelled '1' is now fixed to xi. This new 
system is now defined on the geometry indicated by Figjsj-C (|A| — 1 spins), and governed by the Ising 
Hamiltonian: 

H^'^\(J2 . . .cr|A|) = Hk-l{y^l(J2 . • .cr|A|), 

and its Boltzmann weights obey 



g-^if(2)(a2...an) 



(2) / I A (|A|)/ I X 

= 4_l(^T2|xi) . . . (^|A| |X1 . . . a|A|-l). 



If we now measure the magnetisation at corner '2' for this new system, we get 

m;_i(2|xi) ~ mfc_i(2|xi) = 4-i(t W) " |xi 



^To obtain the last equality, one observes that if an event achieves the maximum, so does the complementary event SJ. 
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The magnetisation m'^_^(2|xi) allows to construct 

Again, this known probability distribution is simple enough that it is possible to draw exactly a sample X2 
according to it. Repeating this reasoning, running along the lattice in the order indicated by the cartoon on 
Figj8]-D, we see that the ability to estimate corner magnetisations combined with Bayes' theorem allows to 
draw sequentially [^according to 

. . . a|A|) = ^fcViKj^A.-ii^ski) . . . (^T|A| |ai . . . a|A|_i). 

Finally, we observe that 



|1- |mfc_i(/|ai...a^_i)| | 

So the condition ( |104| ) leads simply to a condition on the relative precision over the magnetisation. 

Sunmiarising, /6>r any e > 0, temperature /3 and magnetic field h, it is possible to provide an estimate 
Z(h)for the ferromagnetic Ising partition function Z{h) satisfying 

Prob[\Z{h) - Z{h)\ < e Z{h)] > 3/4, (107) 

in a time that scales at most polynomially with (3^ e~^, \h\, and the size of the system if we are able to per- 
form corner magnetisation measurements on related non-homogeneous Ising systems. The required relative 
precision need not be lower than the inverse of some polynomial in /3, and the size of the system. 

F Proof of Theorem Q 

Our starting point is the following direct consequence of the adiabatic theorem, as stated in ll29ll . 



Lemma F.l Let 7 = J^ii^te[0:T]8^P H{f), where gap H(t) denotes the difference between the two lowest 
eigenvalues ofH{t), and let \ ^') denote the quantum state obtained by the continuous evolution induced on 



|$o) by the Hamiltonian family (47). Let also |A| and |£^(A)| denote respectively the number of sites and 
edges of the lattice A. The distance between \^') and the true ground state \G) is at most 5 whenever T 
satisfies 

r>r.A*)^^<^^^i™l^. (108, 

Proof: Let us introduce the parameter s ^t/T. Theorem 2.1 of Ref. |[29ll provides the following sufficient 
condition for adiabaticit)]^ 

in5 w^fjw W-^—frW 

T > n{H,6) = ^ max max{^^^, " ' ' ^^^^ } (109) 
^ 52 o<s<i 7 7 



The order we have chosen has no particular meaning. The reasoning is of course vahd for any labeUing of the sites of the 
lattices. 

^In the following 1 1 A| |oo will denote the operator norm of an operator A, i.e. 1 1^| |oo = sup 3, ^ \fx\\^^ ■ 
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valid for any time-dependent hamiltonian H{t). Adapting this condition to the special case of Hamiltonians 
( |47| ), we see the r.h.s of ( |108| ) certainly upper bounds the r.h.s of ( |109| ). □ 



We wish to discretise the time evolution of our system. Instead of considering the time-dependent evo- 
lution associated with the Hamiltonians H{t), we will deal with L consecutive constant unitary operators, 
Uk = Exp( — i T Hq — i T Hi{kT)) , = . . . L — 1, where we define the discretisation step as 

r = T/L. (110) 

We wish to work with the state l^"*") = Ul-i • • • Z^o|+?'^') rather than with the state |$'). Of course when 
L grows large we expect this substitution to have negligible effect. But we need to be precise and quantify 
the induced error. The following lemma addresses this issue. 

Lemma F.2 The distance between \G) and \ ^~^) is bounded as 



whenever T > T^{H,5). 

Proof: The triangular inequality yields 

II I**) - |G)|| < II 1$') - |G)|| + II 1$*) - |$')ll- (112) 

The first term of the r.h.s of this expression is of course bounded by 5. To bound the second, we use 
Lemma 1 of 1*451, which states that if two time-dependent Hamiltonians Ha{t)^ Hi){t)^0 < t < T differ 
at most by e in operator norm for every t, then the difference between the unitary evolutions they induce, 
Ua{T)Mh{T) satisfy \ \Ua{T) -Ub{T)\\oo < \/2T^. For every t G [0, T], let G {0,...,L- 1} such 
that k{t)T < t < {k(t) + 1)t. Clearly, \\H{t) - H{k{t)r)\\ < r{\h\ • |A| + | J| • \E{A)\). Identifying the 
r.h.s. of this inequality with e and bearing in mind the definition of r, one bounds the second term of the 
r.h.s. of ( |112| ) in the desired way. □ 



Next, we split each unitary Uk into a part that depends only on and a part that depends only on 
Hi{kT)\ for r small enough, each unitary Uk can be safely replaced by the operator 

Uj^^e-'^Ho^-irH,{kr)^ (113) 

Indeed, the Baker-Campbell-Hausdorff identity ll45l implies that 

\\Uk - UkWoo < K{\h\ . |A| + I J| . \E{A)\) . {\h^\ . |A|)T^ (114) 
for some constant K. Then we arrive at the following: 

Lemma F.3 The quantity by which the state Ul-iUl-2 • • • C^o|+x'^') deviates from the true ground state 
ofH^ is at most 



A + +|J|||-|^(A)|) ^ ^^^1^1 ^ ^ i^^^i ^ i^^^^l^ ^ l^^l _ l^l^2_ ^^^^^ 



Proof: The result follows by combining the inequality in Eq. |114| with the Lemmata |F.l|F.2[ □ 
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G Approximation of fidelity overlaps 



In this section we describe how to reconstruct fidelity overlap which is proportional to a partition function 
with complex couplings by sampling from partition functions with real couplings. We begin by rewriting 
Eq. |60] using more compact notation: 

Z(/3)=S(/3) E E E E E E S.i,.2,.3,.4,.5,.ae^-^^'*, (116) 

gi=-ni g2=-n2 g?>=-m gA=-nA gb=-n^ 96=-nQ 

Where /3 = {^\ /3^ p\ /3^ P^} ^ /3, 

m = 712 = 2L|A|, 713 = n4 = 2L^|A|, = L{L - 1)(|A| - 1), uq = L^L' - 1)(|A| - 1). 

and where B0) = Yl^=i ^jW"^) with Bj{P^) = e^^^^^ and c is just a relabeling of c with each index 
ranging from [0, nj] rather than [— nj/2, nj/2]. 
Let us define the polynomial: 

ni 722 ns 724 715 ^6 

^(x) = ^(xi, X2, X3, X4, X5, Xg) = ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ ^ Qi ,22 ,^3 ,^5 7^6-^1^ -^2^ '^3^ '^5 '^6 ' (117) 

Zl=0 i2=0 13=0 24=0 25=0 26=0 

where x = {xi, X2, X3, X4, X5, xe} G M^. Introducing the notation: 

x,-(.) = e-(-) , 

one has the trivial relation: 

pixi{/3+),X2i/3-),X30'+),x4/3'_),X5{f3),xeif3')) = B-\P)Z0) . (118) 

Note that, for physical temperatures, the domain of the polynomial is such that xj > and | |xj 1 1 < 1 for 
j = 1, . . . , 6. We now want to reconstruct the polynomial p{x) from a set of N data values p{xi) with 
x^ G r where F is a certain lattice of points in M^. Although several options are available, in our case the 
polynomial is such that the simplest possible option can be used: a rectangular mesh lattice as: 

r = {Xl^i-^=l, . . . , Xi^i-^=ni+l} X • • • X {X6,i6=l, . . . , Xe^iQ=nQ+l} • 

This is justified by the fact that, as we constructed it, the polynomial p{x) has degree at most nj in Xj 
(j = 1, . . . , 6). This means that p{x) actually lies in the product space 11^^ x • • • Il^g, where 11^ indicates 
the space of univariate polynomials of degree at most n. Explicitely, the data values are written as: 

Pi'^l) = ^(-^1, 21 5 ^2,22 5 ^3,23 5 ^4,24 5 ^5,25 5 ^6,26) =^212223242526 * 

The reconstructed polynomial can then be written as: 

p{x) = P212223242526 ^212223242526 (^) 5 (119) 
2 

where: 

^21222324^5^6 (^) ^ (^i)^2,22 (^2)^3,23 (^3)^,24 (^4)^5,25 (^5)^6,26 (^e) , 
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with: 



'«w= n 

kj = 1 ^''^^ ^' 

kj / ij 



(120) 



We now suppose to have a device that provides an estimate Z{/3) for the partition function, Z{/3), that 
satisfies 

\Z{I5)-Z{(3)\<S , (121) 

and, from this, we want to see how well we can estimate the previously defined overlaps. Since the overlaps 
depend on the analitically continued partition function we are going to show how to reconstruct it. 

From Eq. |118| and Eq. |119| we can write: 



Z{/3*) = 5(/3>(xi(/3;),X2(/3*),X3(/3a),a^4ra,X5(r),X6(/3'*)) 



(122) 



Now the coefficients Piiisis^isie values of the polynomial evaluated at the lattice points Xi, and we 

can use the real temperatures oracle for the partition function in order to write: 



Z(/3^) = S(/3^)5^5-\/3.)Z(/3.)n^.^^^^ 



(123) 



where represents the lattice T transformed in "f3 coordinates": 

f^i — ' /^-,^2 5 /^+,i3 5 /^-,i4 5 A55 /^ie ~ (124) 

= {-logXi,^,, -logX2,i2, -logX3,^3, -logX4,i4, -logX5,i5, -logX6,i6}. 

We also want to make an explicit choice for this lattice: 



^7 

X- = 



Hi + 1 



with: ij = 1, . . . , rij + 1 



(125) 



which clearly satisfies the properties of the rectangular mesh F we stated before. Explicitly the mapping of 
this lattice in the "temperature domain" reads: 



rij + 1 



(126) 



Now, we can write the final formula for the overlap as a function of the estimation of the partition function 
at real temperatures as: 



/ 



1 


r /l-e^i 


L|A| 


r /l-e'2i 


2W 






[y e'4 + 4J 



B0*)^B-\p.)Z0.)fll,^i^{x^(^P*n) . (127) 



i=i 



We are interested in studying how the variance on this quantity scales. We have: 



^/ - 22(i+i'+l)|A| 



(128) 
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From Eq. |121| we have: 

We now study the term by term the quantities in this expression. 



6 . lel^j-r^ 



ni/2 



1 

(l_e2)(L+L0|A| • 



722/2 



n3/2 



714/2 



Then: 



6 



- lli=l yrij+l 

We now turn to each term (^x^{P*^)j \ for each fixed j: 



nj/2 



ii kj =1 \Xj,i.-Xj^kj\ 

kj / ij 



nfe^=i J (Pj COS Oj- 



ij\(nj+l-ij)\ 



where x^{l3^^) = p^e^^^ as can be deduced by looking at Eqs. 55 59 and: 

nl 11 

Hpj^^j)^ / log[(pjcos^^- - (x- — ^))2 + p^sin^/] < --. 



(129) 



(130) 



(131) 



(132) 



(133) 



The last inequality holds for the cases considered by Eq. 59 for nj > 10. Note that 6i — —62 — —^3 = 
O4 = J and 65, ^6 1 since Jt/L, Jt'/L' <C 1 by assumption in the adiabatic mapping. By reassembling 
everything together and from Eq. |129| we get: 



(134) 



with: 



n 



i=i 



(l_e2)2(i+L')|A|22(i+i' + l)|A| 



'3 

, "7+1 , 



"^(„,+l)2"ie-("i+l)/4 



(tj!(n. 



+ l-t,)! 



(l_e2)2(L+L')|A|22(L+L' + l)|A| 



(135) 
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In arriving at the inequality above we used the fact that {pj cos 9j — ^^-^)^ + sin^ 9j > sin^ | for 
j = 1,2,3,4 and {pj cos 9j — ^^;^)^ + pj sin^ 9j > 9j for j = 5, 6, supposing that ^5, ^ as is the 



case. In the "temperature domain" this formula reads: 



n;=i 



/- I i\2n.- -(n.-+l)/4 



r2((n,-+l)e ^^■''i+l)r2((n,-+l)-(n,-+l)e ^^'b' +1) 



'''^■^ gij^8 I ^2 _ ^2)2(L+L^+1)|A|22(L+L^+1)|A| 



(136) 



From this we can get the following condition for the variance on the overelap to be polynomially bounded 
in the system size expressed for generic temperatures: 

S0> < -ji. . (137, 

Explicitely we have: 

vr yj e—T{{nj + l)e-^J + l)T{{nj + 1) - {nj + l)e-^J + 1 

where we used 1 — > ^. By using Stirling approximation we then get: 

\og50)<Y,nA9{Pj) + ^\+Y.\og{nj + l) + Y,9{lij) + K , (139) 

7 = 1 ^ ^ 7 = 1 7 = 1 



803) < e,e, sin^ ^ n ^ ^^""^ ^ ^ ^ ^ » ^ ~ ^ ^ ^ ^ > (138) 



(140) 



where: 

g{/3j) = (l-e-ft)log(l-e-^0-/3ie-^^ -i 
K = 41ogsinf + log6i56i6 . 

This result is telling us how much error we can tolerate in the sampling of the classical partition function 
in order to be able to reconstruct certain quantum overlaps with a precision that scales polynomially in the 
system size. All error values satisfying Eq. 139 allow for such a reconstruction. For this reason, if we want 
to obtain a weaker but more compact result we can chose to state a smaller threshold. We can do this by 
substituting the functions appearing in Eq. [139] with their minimum (see Fig. |G]) getting to: 

log 50) <J2f^ + a)nj + J2 log (nj + !) + &, (141) 

i=i ^ ^ i=i 

where: 

a < min5(/3) ~ -1.6 (142) 
b = K-W . (143) 
In the thermodyamic limit the above formula can be further approximated by: 



(144) 



log S0) < E5=i f I + a) Hj + E;=i log(n,- + 1) + log 9,9^ 
< E?=i(| + «)^i + E5=ilog(n,) + loge5^6 , 
SO that: 

60) < 2\JT){J'T')L{L - 1)L'{L' - l)\AfYley^~^-^r\ (145) 
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I 

QQ, I <N 



I <N 
+ 



Figure 9: Function g[j3) (blue) and + | (red) and f — a with a = min^ g{j3) = —1.6 as a function of 
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